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ABSTRACT 

Every  forecast  should  include  an  estimate  of  its  likely  accuracy,  as  a  measure  of  predictability.  A  new  measure, 
the  first  passage  time  (FPT),  which  is  defined  as  the  time  period  when  the  model  error  first  exceeds  a  predetermined 
criterion  (i.e.,  the  tolerance  level),  is  proposed  here  to  estimate  model  predictability.  A  theoretical  framework 
is  developed  to  determine  the  mean  and  variance  of  FPT.  The  low-order  Lorenz  atmospheric  model  is  taken  as 
an  example  to  show  the  robustness  of  using  FPT  as  a  quantitative  measure  for  prediction  skill.  Both  linear  and 
nonlinear  perspectives  of  forecast  errors  are  analytically  investigated  using  the  self-consistent  Nicolis  model. 
The  mean  and  variance  of  FPT  largely  depends  on  the  ratio  between  twice  the  maximum  Lyapunov  exponent 
(&)  and  the  intensity  of  attractor  fluctuations  (q2),  A  =  2irlq2.  Two  types  of  predictability  are  found:  A  >  1 
referring  to  low  predictability  and  A  <  1  referring  to  high  predictability.  The  mean  and  variance  of  FPT  can  be 
represented  by  the  e-folding  timescales  in  the  low-predictability  range,  but  not  in  the  high-predictability  range. 
The  transition  between  the  two  predictability  ranges  is  caused  by  the  variability  of  the  attractor  characteristics 
along  the  reference  trajectory. 


1.  Introduction 

The  following  is  a  commonly  asked  practical  ques¬ 
tion.  How  long  is  an  atmospheric  (or  oceanic)  model 
valid,  once  it  is  integrated  from  its  initial  state?  To  an¬ 
swer  this  question,  uncertainty  in  atmospheric  models 
must  be  investigated.  It  is  widely  recognized  that  un¬ 
certainty  can  be  traced  back  to  three  factors  (Lorenz 
1984a):  (a)  measurement  errors,  (b)  model  errors  such 
as  discretization  and  uncertain  model  parameters,  and 
(c)  chaotic  dynamics.  Measurement  errors  cause  uncer¬ 
tainty  in  initial  and/or  boundary  conditions.  Discreti¬ 
zation  not  only  causes  truncation  errors  but  also  small- 
scale  “subgrid”  processes  to  be  either  discarded  or  pa¬ 
rameterized.  The  chaotic  dynamics  may  trigger  a  sub¬ 
sequent  amplification  of  small  errors  through  a  complex 
response. 

Traditionally,  the  small  amplitude  stability  analysis 
(linear  error  dynamics)  is  used  with  the  instantaneous 
error  (IE)  growth  rate  and  the  corresponding  e-folding 
timescale  as  the  measures  for  evaluating  the  model  pre¬ 
dictability  due  to  initial  condition  error  (the  first  kind 
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of  predictability).  The  IE  growth  rate  is  usually  esti¬ 
mated  by  either  the  leading  (largest)  Lyapunov  exponent 
or  the  amplification  factors  calculated  from  the  leading 
singular  vectors  (SVs;  e.g.,  Lorenz  1984b;  Dalcher  and 
Kalnay  1987;  Farrell  and  Ioannou  1996a,b;  Vannitsem 
and  Nicolis  1997). 

However,  the  errors  may  grow  to  a  finite  amplitude 
such  as  in  the  medium-range  prediction  (Vukicevic 
1991;  Barkmeijer  1996),  in  the  forecast  using  “imper¬ 
fect”  models  (Bofetta  et  al.  1998;  Palmer  2001),  and  in 
the  models  with  open  boundaries  (Gustafsson  et  al. 
1998;  Chu  1999;  Jiang  and  Malanotte-Rizzoli  1999). 
Thus,  the  linear  assumption  is  no  longer  applicable  and 
the  nonlinear  effect  should  be  considered. 

Both  linear  and  nonlinear  perspectives  of  forecast  er¬ 
rors  can  be  studied  using  a  probabilistic  approach  based 
on  the  analysis  of  the  probability  density  function  (PDF) 
of  the  error  (e.g.,  Benzi  and  Camavale  1989;  Benzi  et 
al.  1999;  Nicolis  1992;  Ehrendorfer  1994a, b;  Molteni 
and  Corti  1998).  The  given  error  PDF  leads  to  the  full 
statistical  description  of  the  error  dynamics.  In  general, 
PDF  satisfies  the  Liouville  equation  which  can  be  solved 
only  by  numerical  methods  even  in  one-dimensional 
cases  (Ehrendorfer  1994a,b;  Nicolis  1992).  A  practical 
approach  to  estimate  the  PDF  of  forecast  error  is  to 
calculate  it  through  ensemble  prediction  (Leith  1974; 
Toth  and  Kalnay  1997). 
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To  simplify  the  probabilistic  analysis  of  forecast  error 
dynamics,  Ivanov  et  al.  (1994)  and  Ivanov  and  Mar- 
golina  (1999)  applied  the  first  passage  time  (FPT;  Gar¬ 
diner  1985)  as  a  measure  of  quantitative  prediction  skill. 
Here,  FPT  is  defined  as  the  time  period  when  the  model 
forecast  error  exceeds  a  predetermined  criterion  s  (i.e., 
the  tolerance  level)  for  the  first  time.  In  this  paper,  we 
plan  to  show  the  advantage  of  using  FPT  in  analyzing 
the  forecast  skill  of  an  atmospheric  model;  to  illustrate 
the  capability  of  such  an  analysis  using  a  low-order 
atmospheric  model;  and  to  improve  the  knowledge  of 
how  the  strength  of  the  initial  error  and  the  variations 
of  model  parameters  affect  the  model  prediction  skill. 

2.  Measure  of  predictability 

a.  Dynamic  constraint 

From  a  low-order  model  to  a  certain  level  of  descrip¬ 
tion  such  as  a  global  atmospheric  spectral  model,  an  N- 
dimensional  vector  x(f)  =  [x(1)(t),  x(2)(f), . . . ,  x(W)(f)]  is 
used  to  represent  the  full  set  of  variables  characterizing 
the  dynamics  of  the  atmosphere.  Let  the  dynamic  law 
be  given  by: 

^=f(x,  t),  (2.1) 

at 

where  f  is  a  functional  including  nonlinear  interactions 
between  different  scales,  a  parameterization  of  subgrid 
scales  and  forcing.  The  aim  of  individual  (oceanic  or 
atmospheric)  prediction  is  to  find  the  solution  of  (2.1) 
with  an  initial  condition 

x(t0)  =  x0. 

Uncertainty  in  atmospheric  (or  oceanic)  models  is 
caused  by  measurement  errors  (initial  and/or  boundary 
condition  errors),  model  discretization,  and  uncertain 
model  parameters.  This  leads  to  the  addition  of  sto¬ 
chastic  forcing  into  atmospheric  (oceanic)  models. 
Moreover,  stochastic  forcing  may  also  be  added  to  the 
forecast  model  in  the  parameterization  of  nonlinear  in¬ 
teractions  (Penland  and  Sardeshmukn  1995;  Leith  1996; 
Palmer  2001). 

For  simplicity,  the  stochastic  forcing  (f ')  is  assumed 
to  be  multiplicative  or  additive  white  noise,  and  (2.1) 
becomes 

dx 

—  =  f(x,  t)  +  f'(x,  t), 
dt 

f '(x,  t)  =  k(x,  t)g(f),  (2.2) 

where  k(x,  t)  and  g (t)  are  the  forcing  covariance  matrix 
{ k(/ }  (dimension  of  IV  X  IV)  and  the  vector  delta-cor¬ 
related  process  (dimension  of  N ),  respectively. 

b.  FPT  for  forecast  error 

Let  x(t)  be  the  reference  solution  which  satisfies  (2.2) 
with  the  initial  condition  x(f0).  The  forecast  error  z  is 
determined  as 


*(0  x'W 


Fig.  1.  Phase  space  trajectories  of  the  reference  solution  x  (solid 
curve),  one  of  the  ensemble  realizations  x',  and  ellipsoids  Sc(t)  and 
S,(r)  centered  at  x.  The  positions  of  x  and  x'  at  time  instance  are 
denoted  by  and  “o”,  respectively.  FPT  is  represented  by  a  time 
at  which  perturbed  trajectory  x'  leaves  the  ellipsoid  Se(t)  at  the  first 
time. 

Z(t)  =  x(f)  -  x'(f), 

where  x'(f)  is  an  individual  prediction  corresponding  to 
perturbing  initial  condition  and/or  forcing  f  +  f '.  The 
rmse  is  used  to  estimate  the  model  forecast  skill: 

J  =  (zTAz),  (2.3) 

where  A  is  the  weight  matrix,  the  superscript  T  denotes 
the  transpose  operator,  the  bracket  represents  the  en¬ 
semble  average  over  realizations  generated  by  stochastic 
forcing,  uncertain  initial  conditions,  and  uncertain  mod¬ 
el  parameters  in  (2.2). 

To  quantify  FPT,  we  first  define  two  model  error  lim¬ 
its.  First,  the  forecast  error  cannot  be  less  than  a  min¬ 
imum-scale  8,  which  depends  on  the  level  of  thermal 
and  intrinsic  noises  existing  in  the  model.  In  general, 
this  scale  functionally  depends  on  the  square  root  of  the 
subgrid-scale  kinetic  energy  Em.  The  simplest  param¬ 
eterization  is  8  proportional  to  Em,  with  the  propor¬ 
tionality  depending  on  the  selection  of  the  subgrid-scale 
parameterization  (Aurell  et  al.  1996).  Second,  the  fore¬ 
cast  error  cannot  be  more  than  a  maximum-scale  (tol¬ 
erance  level)g. 

The  ratio  between  the  maximum  and  minimum  limits 
is  large 

z,  =  e/8  »  1,  (2.4) 

in  atmospheric  circulation  models.  Two  ellipsoids  Ss  and 
Ss  are  defined  by  J  =  e 2  and  J  =  82,  with  the  center 
at  x(t)  (Fig.  1).  The  model  rmse  ( J )  is  bounded  by  the 
two  limits 

82  <  J  <  s2.  (2.5) 

FPT  is  a  stochastic  time  period  rfor  which  x'(f)  leaves 
a  domain  limited  by  Ss  for  the  first  time.  From  physical 
point  of  view,  it  determines  the  low  limit  of  the  pre- 
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dictability  time  because  we  neglect  the  impact  of  res¬ 
toration  of  model  prediction  skill  a  priori,  which  occurs 
in  large  atmospheric  models  (Anderson  and  van  der 
Dool  1994). 

To  calculate  the  moments  of  FPT,  the  probability  den¬ 
sity  of  prediction  (PDP)  is  P(t0,  z0,  t  -  f0),  which  is 
the  probability  of  the  time  period  t  —  t0  that  the  forecast 
error  is  less  than  e  when  the  initial  error  is  z0  (Ivanov 
and  Margolina  1999).  Thus,  the  stability  analysis  on  the 
reference  solution  x  (Fig.  1)  using  FPT  should  be  clas¬ 
sified  as  the  probabilistic  (stochastic)  stability  analysis 
(Has’minskii  1980). 

Use  of  FPT  and  PDP  has  the  following  advantages 
in  evaluating  model  predictability  skill.  First,  FPT  is  a 
measure  of  model  forecast  skill  for  both  linear  and  non¬ 
linear  error  regimes.  Second,  the  moments  of  FPT  char¬ 
acterize  high-order  variability  of  the  prediction  skill. 
This  variability  is  determined  only  from  the  knowledge 
of  a  single-point  distribution  function,  PDP.  Third,  PDP 
and  naturally,  the  moments  of  FPT,  can  be  effectively 
calculated  by  ensemble  prediction  techniques.  Such  a 
calculation  consumes  less  computer  resources  because 
PDP  is  determined  only  inside  the  ellipsoid  Se,  rather 
than  in  the  whole  phase  space.  Fourth,  an  iterative  tech¬ 
nique  (Ivanov  et  al.  1999)  and  analytical  approach  (Chu 
et  al.  2002)  are  developed  to  find  the  structure  of  PDP. 
Analytical  estimate  of  PDP  is  possible  even  for  multi¬ 
dimensional  error  models. 

3.  Pontryagin-Kolmogorov  equation 

The  temporal  evolution  of  PDP,  for  the  model  (2.2), 
satisfies  the  Pontryagin-Kolmogorov  equation  (PKE; 
Pontryagin  et  al.  1969)  with  the  first-order  partial  dif¬ 
ferentiation  with  respect  to  time,  and  the  second-order 
partial  differentiations  with  respect  to  initial  condition 
errors, 

~-(f-  -hk^-  =  0 

dt  \L  dt  )dz°  2  “k,Jdz°z° 

i,  j,  Z  =  1,  ....  A,  (3.1) 

where  the  coefficients  ku  are  the  components  of  the  forc¬ 
ing  covariance  matrix  k(x,  t)  and  (zf,  z°,  •  •  •  z°)  are  the 
components  of  the  initial  error  z0. 

One  initial  condition  and  two  boundary  conditions 
(with  respect  to  z0)  are  needed  to  solve  (3.1).  The  initial 
condition  of  (3.1)  is  given  by 


P(t0,  z0,  t  -  t0)  =  0  at  J  =  e2,  (3.3) 

which  is  called  the  absorbing  type  boundary  condition. 
It  is  clear  that  the  forecast  error  can  not  be  less  than  8, 
which  leads  to 

kuTi  =  0  at  J  =  8\  (3.4) 

dzj 

which  is  called  the  reflecting  boundary  condition.  In 
general,  the  reflecting  and  absorbing  boundary  condi¬ 
tions  can  be  simultaneously  applied  on  Sc  and  Ss. 

The  kth  moment  (k  =  1,  2,  . . .)  of  FPT  is  calculated 
using  PDP  by 

Tt(z0)  =  k  I  P{t0,  z0,  t  -  t0)(t  -  toy-'  dt, 

J  to 

k  =  1,  . . . ,  oo.  (3.5) 

The  mean  and  variance  of  FPT  can  be  calculated  from 
the  first  two  moments: 

<t>  =  r,.  (3.6) 

(fir2)  =  r2  -  r2,  (3.7) 

where  the  bracket  denotes  the  ensemble  average  over 
realizations  generated  by  stochastic  forcing.  The  mean 
and  variance  of  FPT  are  both  local  characteristics  of  the 
prediction  skill  since  they  depend  on  position  z0  on  the 
reference  trajectory  x. 

The  prediction  skill  of  atmospheric  models  is  cur¬ 
rently  verified  using  the  perfect  model  concept  with 
uncertainty  appearing  only  in  the  initial  conditions. 
Such  an  uncertainty  is  easily  taken  into  account  in 
our  approach  through  additional  averaging  of  (3.6) 
and  (3.7)  with  respect  to  an  ensemble  of  initial  per¬ 
turbations. 

4.  A  low-order  atmospheric  model 

Determination  of  prediction  skill  of  a  low-order  at¬ 
mospheric  model  proposed  by  Lorenz  (1984a)  is  taken 
as  an  example  to  demonstrate  the  usefulness  of  our  ap¬ 
proach.  The  Lorenz  system  is  the  simplest  possible  mod¬ 
el  capable  of  representing  an  unmodified  or  modified 
Hadley  circulation,  determining  its  stability  and,  if  it  is 
unstable,  representing  a  stationary  or  migratory  distur¬ 
bance.  The  model  consists  of  the  three  ordinary  non- 
dimensional  differential  equations: 


P(t0,  z0,  0)  =  1,  (3.2) 

which  indicates  that  the  vector  of  initial  error  z0  is  al¬ 
ways  less  than  the  given  tolerance. 

The  model  error  vector  z  is  bounded  by  two  limits: 
minimum  error  limit  (noise  level  8)  and  maximum  error 
limit  (tolerance  level  e).  If  the  initial  error  vector  z0 
reaches  the  tolerance  level  [i.e.,  z0  hits  the  boundary  of 
Sc  (/„)],  the  model  loses  prediction  capability  at  the  first 
place: 


— x\  —  x\  —  axl  +  aF, 
x,x2  —  bx,x3  —  x2  +  G, 

bx,x2  +  xtx3  —  x3,  (4.1) 

Jf;  .Ujf  ntMtfbfl'  •  i 

where  t  is  the  time  scaled  by  five  days;  jc,  represents 


dxx 

dt 

dx2 

~dt 

dXy 

it 
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the  intensity  of  the  symmetric  globe-encircling  westerly 
wind,  and  also  the  poleward  temperature  gradient  which 
is  assumed  to  be  in  permanent  equilibrium  with  the 
wind;  x2  and  x3  are  the  cosine  and  sine  phases,  respec¬ 
tively,  of  a  chain  of  superposed  large-scale  eddies, 
which  transport  heat  at  a  rate  proportional  to  the  square 
of  their  amplitude,  and  transport  no  angular  momentum 
at  all.  The  terms  aF  and  G  represent  symmetric  and 
asymmetric  thermal  forcing.  The  terms  xtx2  and  xtx3 
describe  amplification  of  eddies  through  interaction  with 
the  westerly  currents.  The  displacement  of  eddies  by  the 
westerly  wind  is  parameterized  using  the  terms  —bx1x3 
and  bxIx2.  The  parameter  a  is  the  damping  coefficient. 

The  traditional  stability  analysis  on  the  Lorenz  model 
(4.1)  leads  to  three  Lyapunov  exponents.  Among  them, 
only  one  is  positive.  Projection  of  the  three-dimensional 
forecast  error  vector  onto  the  unstable  manifold  leads 
to  a  self-consistent  model  (Nicolis  1992) 

^  =  (o'  -  g£)£  +  v(t)£,  fU*  =  £o> 

€  6  [0,  »).  (4.2) 

where  £  is  the  nondimensional  amplitude  of  error,  and 
g  is  a  nonnegative,  generally  time-independent  nonlin¬ 
ear  parameter  whose  properties  depend  on  the  under¬ 
lying  attractor.  Equation  (4.2)  is  written  in  Ito  form.  The 
tangent  approximation  of  error  growth  leads  to  g  =  0. 
The  eigenvector  £  is  associated  with  the  positive  Lya¬ 
punov  exponent  (cr).  The  term  v(t)£  is  a  specially  chosen 
stochastic  forcing  with  zero  mean  and  pulse-type  var¬ 
iance 

<v«>  =  0,  =  q2A(t  -  t'),  (4.3) 

where  the  bracket  is  defined  as  the  ensemble  mean  over 
realizations  generated  by  the  stochastic  forcing,  A  is  the 
Delta  function,  and  q2  is  the  intensity  of  attractor  fluc¬ 
tuations  modeled  as  multiplicative  noise  forcing. 

The  combination  of  the  free  model  parameters  cr,  g, 
and  q2  affects  the  model  prediction  skill.  Nicolis  (1992) 
used  the  following  values: 

a  =  0.64,  g  =  0.3,  q2  =  0.2.  (4.4) 

In  reality  the  dynamical  characteristics  of  the  atmo¬ 
sphere  vary  considerably  between  different  synoptic  sit¬ 
uations.  Model  forecast  skill  depends  on  various  factors 
such  as  season,  location,  and  boundary  conditions.  Ev¬ 
idence  shows  variability  of  forecast  skill  in  operational 
models  such  as  the  European  Centre  for  Medium-Range 
Weather  Forecasts  (ECMWF)  model.  Therefore,  the  pa¬ 
rameters  cr,  g,  and  q2  should  be  time  and  attractor  de¬ 
pendent.  This  problem  is  not  addressed  here.  However, 
in  contrast  to  Nicolis  (1992)  we  assume 

0.2  <  < r  <  0.64,  g  =  0.3, 

0.01  <  q2  <  0.6  (4.5) 

in  this  study.  Since  the  Nicolis  model  is  analytical,  the 
parameter  8  should  be  small  and  z,  =  s/8  =  105  is 


assumed.  This  should  be  satisfactory  for  the  quality 
analysis  of  the  prediction  skill. 

Use  of  the  self-consistent  model  (4.2)  has  several  ex¬ 
plicit  advantages.  First,  it  allows  analytical  study  on 
linear  and  nonlinear  perspective  of  forecast  error.  Sec¬ 
ond,  the  methodology  developed  by  Nicolis  (1992)  can 
be  generalized  to  more  realistic  atmospheric  models. 
Third,  the  model  contains  several  dynamical  regimes  of 
forecast  error  behavior.  Obviously,  their  analytical  study 
may  be  useful  for  the  interpretation  of  results  obtained 
by  large  atmospheric  models. 


5.  Mean  and  variance  of  FPT 

For  the  self-consistent  error  model  (4.2),  the  mean 
and  variance  of  FPT  can  be  computed  without  the  de¬ 
termination  of  PDP.  PDP  for  the  Nicolis  model  (4.2) 
satisfies  the  following  PKE  (Pontryagin  et  al.  1969;  Gar¬ 
diner  1985;  Ivanov  et  al.  1994): 


dP  dP  1  d2P 

~  ~  (<r  -  gZo)ZoJT  -  =  0,  (5.1) 

ot  og0  L 


d£l 


with  the  initial  condition 

P(to,  £»  0)  =  1, 

and  boundary  conditions 

P(t0,  £0,  t  —  f0)lfo=e  — 
dP(t0,  £0,  t  -  f0)| 


(5.2) 


d£0 


=  0. 


(5.3) 


£o=S 


Multiplying  (5.1)  by  (t  —  t0)k  and  then  integrating  with 
time  from  f0  to  00 ,  we  obtain  a  set  of  equations  for  the 
moments  r k  (k  —  1,  2,  . . .): 

«*.  -  tufe  +  tf%  - 

6  <  <  e.  (5.4) 

The  boundary  conditions  for  r*  also  are  obtained 
through  integrating  (5.3)  with  time 


=  0, 


dTk 
d£ o 


=  0. 


(5.5) 


The  mean  FPT  is  the  same  as  the  first  moment  r,, 
and  the  variance  of  FPT  is  computed  from  the  first  and 
second  moments  [see  (3.7)].  Simple  mathematical  ma¬ 
nipulations  give  the  following  equations: 

,  ,  „,d{T)  ,  q2£l  d2{r)  ,  J  ^ 

(ofo  -  8to)-jr  +  ~T"“777  =  _1  and  (5-6) 
d£ o  2  d£l  t. 

„  d(8T2)  q2£20d2(8T2)  „ Jd(r)Y 

(»f”  -  +  —  -« • 

(5.7) 


2864 


JOURNAL  OF  THE  ATMOSPHERIC  SCIENCES 


Volume  59 


iii  io  vcitmA  jam!: 


Fig.  2.  Comparison  of  prediction  skill  between  iterative  (with  different  iteration  parameter  Aja)  and 
numerical  solutions  for  (a)  mean  FPT  with  A  =  6.4,  (b)  variance  of  FPT  with  A  =  6.4,  (c)  mean  FPT  with 
A  =  1,  and  (d)  variance  of  FPT  with  A  =  1. 


6.  Iterative  solutions 

An  iterative  method  (Ivanov  et  al.  1999),  whose  math¬ 
ematical  detail  is  listed  in  appendix  A,  is  proposed  to 
solve  (5.6)  and  (5.7).  The  iteration  parameter  (A /a)  is 
defined  as  the  inverse  of  the  number  of  iteration  (7); 
that  is,  A/u,  =  1/7. 

It  is  clear  that  the  mean  and  variance  of  FPT,  as  well 
as  the  iteration  parameter  depend  on  the  model  param¬ 
eters  (o-  and  q)  and  the  tolerance  level.  To  illustrate  how 
fast  the  convergence  of  the  iterative  process  is,  the  it¬ 
erative  solutions  ((t),  (St2))  are  compared  with  the  nu¬ 
merical  solutions  ((r)num,  <ST2)num),  obtained  using  the 
finite  difference  method  with  the  second-order  accuracy 
(Press  et  al.  1992).  Two  sets  of  calculations  are  con¬ 
ducted  (Fig.  2)  with  the  model  parameters,  (a)  In  case 


1,  g  =  0.3,  tr  =  0.64,  q2  =  0.4  (A  =  3.2),  and  e  =  0.1. 
(b)  In  case  2,  g  =  0.3,  a  =  0.2,  q2  =  0.4  (A  =  1),  and 
e  =  0.5.  Here,  A  =  2 <rlq2. 

The  root-mean-square  errors  between  iterative  and 
numerical  solutions. 


((t)  -  (t)„, 
7),  =  rmse| 


<T>„, 


((St2)  -  (St2)\ 


(6.1) 

(6.2) 


are  used  to  represent  the  accuracy  of  the  iteration  pro¬ 
cedure.  For  case  1,  77,  =  0.85%  and  D2  =  5.6%  with 
the  iteration  parameter  A /la  =  1/2;  and  D,  =  0.03%  and 
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D2  =  0.25%  with  the  iteration  parameter  A/x  =  1/3  (Figs. 
2a, b).  For  case  2,  D,  =  1.7%  and  D2  =  5.2%  with  the 
iteration  parameter  A/x  =  1/2;  and  D,  =  0.8%  and  D2 
=  2.2%  with  the  iteration  parameter  A/x  =  1/3  (Figs. 
2c,d).  Thus,  with  the  iteration  parameter  A/x  =  1/3,  the 
iterative  solutions  are  quite  close  to  the  numerical  so¬ 
lutions,  and  therefore,  can  be  used  for  evaluating  the 
model  prediction  skill. 

7.  Characteristics  of  prediction  skill 

The  model  prediction  skill  (5)  that  is  represented 
by  the  mean  and  variance  of  FPT  strongly  depends 
on  the  model  parameters  such  as  A  in  the  Nicolis  mod¬ 
el.  If  a  critical  value  A0  is  found  where  dS/dA  is  dis¬ 
continuous,  S  may  behave  differently  between  A  > 
A0  and  A  <  A0.  Low  (high)  S  values  in  one  side  of 
A0  is  referred  to  low  (high)  predictability.  Low  pre¬ 
dictability  (dynamical  chaos)  may  coexist  with  high 
predictability  (dynamical  behavior)  even  for  the  same 
attractor  (Ziehmann  et  al.  2000).  In  the  Nicolis  sys¬ 
tem,  low  predictability  is  found  for  A  >  1,  and  high 
predictability  is  found  for  A  <  1. 


a.  Low  predictability  (A  >  1) 

Low  predictability  is  characterized  by  low  values  of 
mean  and  variance  of  FPT,  high  sensitivity  to  initial 
condition  errors,  and  usually  occurs  in  chaotic  dynam¬ 
ics.  This  case  has  been  numerically  studied  by  Nicolis 
(1992)  for  A  =  6.4. 

1)  Linear  regime 

For  the  linear  regime  ( g  =  0),  the  mean  and  variance 
of  FPT  have  the  analytical  expressions 


<T)[0] 


Z 


q\ A  -  1) 


ln(z)  - 


A  -  1 


zr(A-l)(2X- 


£>’ 


e 

S' 


1)  , 

(7.1) 


Since  z,  »  1  and  A  >  1,  Eq.  (7.1)  can  be  rewritten  by 

<t>™  -  - . 2-  -2,  ln(z),  (7.2) 

(2cr  -  q 2) 

which  is  the  classical  e-folding  time.  Forecast  accuracy 
weakly  depends  on  8  and,  naturally,  on  the  model  sub- 
grid-scale  parameterizations.  Furthermore,  (t)[01  is 
bounded  if  S  is  nonzero. 

Equation  (7.2)  further  shows  that  (t)[0]  increases  with 
increasing  intensity  of  attractor  fluctuations  q 2 .  This  is 
consistent  with  the  Benzi  and  Camevale  (1989)  result 
obtained  for  the  tangent  model  (4.2)  that  the  most  prob¬ 
able  growth  of  forecasting  error  is  depressed  as  the  at¬ 
tractor  fluctuations  increase. 

Similar  to  (7.1),  the  variance  of  FPT  is  calculated  by 


<5t2>I°1  ~  4n  8  ns  ln<z>’  (7-3> 

q\ A  -  l)3 

which  indicates  that  the  variance  of  FPT  also  has  an  e- 
folding  dependence  on  z.  However,  the  generality  of  this 
analytical  result  (i.e.,  e-folding  variance  of  FPT)  into 
operational  atmospheric  models  should  be  further  in¬ 
vestigated.  The  mean  and  variance  of  FPT  are  calculated 
with  three  different  tolerance  levels  (s  =  0.02,  0.2,  and 
2)  and  various  values  of  initial  error  £0.  Both  mean  and 
variance  of  FPT  reveal  e-folding  type  with  z  when  £0 
<  0.5  (Fig.  3). 


2)  Nonlinear  regime 

For  the  nonlinear  regime  ( g  ¥=  0),  the  mean  FPT  is 
obtained  using  the  iterative  method  illustrated  in  ap¬ 
pendix  A: 

<t)  »  (t)™  + - - - (/?!1)  +  fl'21  +  R 131  +  fl[41), 

3</2(A  -  1)  »?- 

<7.4) 

where  the  functions  /?m,  R[2],  Rl3\  and  /?141  are  deter¬ 
mined  as 


Rm  =  A,£f(l  -  z-1), 

P3£3 

R™  =  A3- — (1  -  z-3),  and 
q6 

g4e4 

RM  =  A4S—  (1  -  z~4). 


g2e2 

R™  =  A,5— (1  -  z-2), 

ZJ4 


with  the  coefficients 
5 


16 


A‘0-  A,  A2  -  A(A  +  iy 


A,  = 


A4  = 


5.9 


,  and 


A(A  +  1)(A  +  2) 

_ 4/7 _ 

A(A  +  1)(A  +  2)(A  +  3)' 


As  an  example,  for  A  =  3.2,  the  coefficients  are: 

A,  =  1.5626,  A2  =  0.678,  A3  =  0.030, 

A4  =  0.002, 

which  indicates  a  fast  convergence  of  the  iteration 
process  (fast  reduction  values  from  A,  to  A 4).  If  A  = 
6.4  [the  parameters  take  on  the  same  values  as  in 
Nicolis  (1992)],  the  convergence  of  the  iteration  is 
even  faster. 

The  expression  for  the  variance  of  FPT  is  quite  cum¬ 
bersome,  and  is  not  listed  here.  Figure  3  shows  the 
comparison  of  the  mean  and  variance  of  FPT  between 
linear  and  nonlinear  regimes  with  three  different  tol¬ 
erance  levels  (e  =  0.02,  0.2,  and  2)  and  A  =  3.2.  For 
small  tolerance  levels  (e  =  0.02,  0.2),  there  is  almost 
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Fig.  3.  Comparison  of  prediction  skill  between  linear  and  nonlinear  regimes  for  (a)  mean  FPT  with  e 
=  0.02,  (b)  variance  of  FPT  with  e  =  0.02,  (c)  mean  FPT  with  e  =  0.2,  (d)  variance  of  FPT  with  e  = 
0.2,  (e)  mean  FPT  with  e  =  2,  and  (f)  variance  of  FPT  with  e  =  2. 


no  difference  between  linear  and  nonlinear  regimes.  The 
difference  increases  with  the  tolerance  level.  As  the  tol¬ 
erance  level  increases  to  s  =  2,  the  difference  becomes 
evident. 

For  a  small  initial  error  (£0  <C  1),  (7.4)  becomes 


<T)  5=3  ^T^(ln(z)  +  Vx8e 

32 


3A(A  +  1)  q 


g2s2 

'  T.  TV* 


7  i' 


,  (7.5) 


which  indicates  that  the  increase  of  the  mean  FPT 
with  increasing  q2\  that  is  to  say,  the  mean  FPT,  is 
larger  for  the  nonlinear  regime  than  for  the  linear 
regime. 


For  a  large  initial  error  £0,  the  dependence  of  the 
high-order  iterations  on  the  stochastic  forcing  (q2)  is 
complicated.  For  example,  the  impact  of  the  third  it¬ 
eration  is  given  by 


U  = 


3<72(A  -  1) 


-RW 


18.8g4e4 


3q'°(k  -  1)A(A  +  1)(A  +  2)(A  +  3) 


(1  -  z-4). 


Three  U  values  (for  q1  =  0.01,  0.02,  and  0.2),  and  in 
turn  the  two  ratios,  are  computed: 

U(q2  =  0.01  )/U(q2  =  0.02)  =  1.04, 

U(q2  =  0.02  )/U(q2  =  0.2)  =  0.53, 
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I  '  ■  F?  i,  _  J ti 1  ■  oi 

Fig.  4.  Dependence  of  the  prediction  skill  on  A  for  (a)  mean  FPT  with  e  =  0.1,  (b)  variance  of  FPT  with 
e  =  0.1,  (c)  mean  FPT  with  e  =  0.5,  and  (d)  variance  of  FPT  with  e  =  0.5. 


which  indicates  that  U  is  reduced  when  q2  increases 
slightly  from  q2  =  0.01  to  q2  =  0.02;  and  then  is  in¬ 
creased  when  q2  increases  drastically  from  q2  =  0.02 
to  q2  =  0.2. 

b.  Transition  to  high  predictability  (\  =  1) 

The  point  A  =  1  is  a  threshold  between  two  kinds  of 
prediction  skill.  Here,  the  mean  and  variance  of  FPT 
cannot  be  described  by  the  e-folding  type  dependence 
on  z  even  for  infinitesimal  and  small  initial  errors.  They 
depend  on  the  noise  level  S.  To  show  a  change  of  the 
model  predictability  at  A  =  1,  the  mean  and  variance 
of  FPT  are  calculated  using  (5.6)  and  (5.7)  with  various 


tolerance  levels  and  initial  errors  for  linear  and  nonlinear 
regimes. 

The  two  different  predictability  ranges  (low  and  high) 
can  be  identified  from  the  dependence  of  mean  and  var¬ 
iance  of  FPT  on  A  value  (Fig.  4).  Notice  that  the  slope 
of  curves  and  functional  dependences  for  the  mean  and 
variance  of  FPT  have  a  drastic  change  at  around  A  = 
1,  that  is,  the  low-predictability  range  (low  mean  and 
large  variance  of  FPT)  to  high-predictability  range  (high 
mean  and  large  variance  of  FPT). 

1)  Linear  regime 

For  the  linear  regime  ( g  =  0),  the  mean  and  variance 
of  FPT  are  given  by 
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Fig.  5.  Dependence  of  prediction  skill  on  f0/e  for  linear  and  nonlinear  regimes  with  A  =  1  (transition 
between  low  and  high  predictability)  for  (a)  mean  FPT  with  e  =  0.01,  (b)  variance  of  FPT  with  e  =  0.01, 
(c)  mean  FPT  with  e  =  0.1,  (d)  variance  of  FPT  with  e  =  0.1,  (e)  mean  FPT  with  e  =  0.5,  and  (f)  variance 
of  FPT  with  e  =  0.5. 


<T>‘°1  =  ~^[(A,)2  -  (A2)2] 

=  ~2  ln(  — ) [ln(ezo)  -  InS2], 

q2  \zo 


(7.6) 


of  FPT  (proportional  to  S2)  decreases  with  decreasing 
82.  Therefore,  the  reduction  of  the  noise  level  only  de¬ 
creases  the  variance  of  FPT  considerably,  but  not  the 
mean  FPT. 


<8t2>'01  =  ^[(Al)4  _  (A2)4]’ 
where 


(7.7) 


A!  =  ln(z,),  A2  =  ln(z2),  z2  =  £0/S. 


Several  interesting  facts  emerge  from  (7.6)  and  Fig.  5: 
(a)  the  mean  FPT  does  not  coincide  with  the  e-folding 
timescale  [see  (7.2)];  (b)  the  mean  FPT  increases  with 
decreasing  8  logarithmically  as  —  ln(S2);  (c)  the  variance 


2)  Nonlinear  regime 

iV  ■>  ■■  .-.vVAtrs  -.,W«  . 

For  nonlinear  regime  (g  ^  0),  the  mean  FPT  is 
obtained  using  the  iterative  method  illustrated  in  ap¬ 
pendix  A: 

vufi-  -pb  aqw  Mitfbioi  --1  agi  yd  iv 

<T>  -  <r><°>  +  +  R™  +  R™  +  RW),  (7.8) 

.  r?  !  ,,  t 

where 
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Rm  =  5i~[(Ai)3  -  (A2)3], 

q2 

r™  t(A,)4  -  (a 2n 

q* 

Rm  =  fi3^[(A1)5  -  (A2)5], 

q 

RW  =  b4^[(a ,)«  -  (a2)«], 
q 8 

and  the  coefficients  B{,  B2,  B3,  and  fi4  are  given  by 
B{  =  2.66;  B2  =  1.33;  B3  =  0.36; 

B4  =  0.04. 

Figure  5  shows  the  comparison  of  the  mean  and  var¬ 
iance  of  FPT  between  linear  and  nonlinear  regimes  with 
three  different  tolerance  levels  (e  =  0.01,  0.1,  and  0.5). 
For  the  small  tolerance  level  (e  =  0.01)  the  nonlinear 
effect  on  both  the  mean  and  variance  of  FPT  is  negli¬ 
gible.  For  0.01  <  e  <  0.1,  the  nonlinear  effect  is  neg¬ 
ligible  on  the  mean  FPT  too  but  not  on  the  variance  of 
FPT.  For  b  =  0.5,  the  prediction  skill  difference  between 
the  tangent  and  nonlinear  error  models  becomes  very 
evident.  For  the  same  initial  error,  the  mean  FPT  is  larger 
for  the  nonlinear  regime  than  for  the  linear  regime, 
which  indicates  the  enhancement  of  the  predictability 
due  to  nonlinearity. 

c.  High  predictability  (\  <  1) 

As  A  <  1,  the  mean  and  variance  of  FPT  increase 
with  decreasing  A  rapidly  (Fig.  4),  which  leads  to  the 
high-predictability  range. 

1)  Linear  regime  (g  =  0) 

For  A  — »  1,  the  mean  and  variance  of  FPT  are  given  by 


<t)m  - 

71 — ^77ia>)1-A[l  -  and 

(7.9) 

(1  -  A )V 

(gT2)lOJ  * 

1 

(7.10) 

Equation  (7.9)  demonstrates  that  the  mean  FPT  does 
not  coincide  with  the  e-folding  timescale  (7.2)  even 
for  infinitesimal  and  small  amplitude  errors.  In  ad¬ 
dition  it  is  noted  that  for  high  predictability  and  z  » 
1,  the  mean  and  variance  of  FPT  have  power  depen¬ 
dence  of  z,  with  scaling  exponents  equaled  to  1  —  A 
and  2(1  -  A),  respectively. 


(t>  =  Mm  + - — - (Rm  +  Rm  +  Rm  +  Rwy 

r-o.fs-s  2>q2(\  -  A) 

(7.11) 

where 

i  A  :  ' 

Rm  =/t—n  -  z2(1~A)),  =  C2—(\  -  Zw-A>), 

/  q2  q3 

Rm  =  c3^(i  -  z*1-*)), 

q 

„3e3 

R™  =  C4^—(l  -  z5«-A>), 

q 5 

and  Cl9  C2 ,  C3,  and  C4  are  the  coefficients  given  by 

Ci  = - -  Co  = - *  c,  — - ■ 

1  2  y  A  (3  -  A)  -!  9(4  -  A) 

7  0 

c  =  32 

4  135(5  -  A)’ 

For  A  =  0.5,  they  are 

C,  =  0.67,  C2  =  1.2,  C3  =  0.13,  C4  =  0.05, 

which  indicates  a  slower  convergence  of  the  iteration 
process  (reduction  values  from  C,  to  C4)  than  in  the 
low-predictability  regime  (A  »  1)  represented  by  the 
reduction  values  from  A,  to  A4. 

Figure  6  shows  the  comparison  of  the  mean  and 
variance  of  FPT  between  the  linear  and  nonlinear  er¬ 
ror  models  for  three  different  tolerance  levels  (e  = 
0.01,  0.1,  and  0.5).  For  the  small  tolerance  level  (s 
=  0.01),  the  nonlinear  effect  on  both  mean  ^nd  var¬ 
iance  of  FPT  is  negligible.  As  the  tolerance  level  in¬ 
creases  (i.e.,  e  >  0.01),  the  mean  FPT  is  larger  in  the 
nonlinear  regime  than  in  the  linear  regime  for  the 
same  initial  condition  error  vector,  which  indicates 
the  enhancement  of  the  predictability  due  to  nonlin¬ 
earity.  Besides,  both  solid  and  dotted  curves  in  Fig. 
6  reveal  a  reverse  T -type  shape,  which  indicates  that 
(t)  and  (St2)  depend  weakly  on  the  value  of  initial 
error  in  both  linear  and  nonlinear  regimes  if  £0/e 
0.01. 

8.  Ensemble  effect  of  initial  uncertainty  on  the 
mean  FPT 

The  error  of  an  ensemble  forecast  is  usually  smaller 
than  the  mean  error  of  individual  forecasts  comprising  the 
ensemble  (Leith  1974;  Murthy  1991).  One  may  ask:  What 
is  the  ensemble  effect  of  initial  uncertainty  with  finite 
amplitudes  on  the  prediction  skill?  The  classical  stability 
analysis  associated  with  Lyapunov  exponents  or  SVs  is 
not  able  to  answer  this  question  since  it  is  in  general  only 
for  small  initial  errors.  The  FPT  approach  can  answer  this 
question,  and  the  low  predictability  range  (A  >  1)  is  taken 
as  the  example  for  illustration.  To  do  so,  let  the  initial 
condition  be  distorted  by  noises  with  unknown  intensity 


2)  Nonlinear  regime  (g  0) 

For  the  nonlinear  regime,  the  mean  and  variance  of 
FPT  are  given  by 
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of  |0.  According  to  Leith  (1974)  additional  noises  with 
zero  mean  and  variance  of  4f 2  are  added  to  f0.  The  mean 
FPT  is  averaged  over  ensemble  realizations  composed  of 
individual  samples  with  noises  homogeneously  distributed 
within  the  interval  [|0  —  f,  |0  +  f], 

a.  Infinitesimal  and  small  initial  errors 

Independent  of  the  tolerance  level,  we  obtain  the  fol¬ 
lowing  expression: 

«r(lo)>)  >  <r(io)>,  (8.1) 

for  infinitesimal  and  small  initial  errors.  Here,  the  single 
bracket  indicates  an  averaging  over  statistics  of  sto¬ 
chastic  forcing  and  the  double  bracket  denotes  an  av¬ 
eraging  operator  over  statistics  of  both  stochastic  forcing 
and  the  ensemble  of  initial  realization.  The  mathematical 


,  except  for  A  =  0.5. 

detail  is  listed  in  appendix  B.  The  inequality  (8.1)  shows 
that  the  ensemble  mean  of  forecasts  will,  on  average, 
have  longer  FPT  (smaller  model  error)  than  the  mean 
error  of  individual  forecast  comprising  the  ensemble. 
This  agrees  with  the  classical  result  (Leith  1974;  Mur¬ 
phy  1991). 

b.  Finite  amplitude  initial  errors 

For  finite  amplitude  initial  errors,  the  ensemble  mean 
FPT  may  be  smaller  than  the  mean  FPT  for  an  individual 
prediction.  Therefore,  if  the  initial  uncertainty  consists 
of  small  and  finite  amplitude  errors  and  if  the  tolerance 
level  is  large,  the  additional  averaging  over  the  ensemble 
of  the  initial  errors  may  deteriorate  the  model  prediction 
skill. 
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9.  Conclusions 

This  study  demonstrates  the  capability  of  the  FPT 
approach  in  evaluating  model  forecast  skill.  Two  kinds 
of  predictability  skill  are  analytically  identified  for  the 
Lorenz  (1984a)  model  based  on  the  parameter,  A  =  2<r/ 
q 2:  low  predictability  if  A  >  1,  and  high  predictability 
if  A  <  1.  Here,  cr  and  q  are  the  maximum  Lyapunov 
exponent  and  the  intensity  of  attractor  fluctuations,  re¬ 
spectively.  The  low-  to  high-predictability  transition  oc¬ 
curs  at  A  =  1 .  The  transit  from  low  to  high  predictability 
is  realized  for  A  =  1.  Both  kinds  of  predictability  skill 
coexist  for  small  initial  error  £0  <C  1. 

The  low  predictability  (A  >  1)  is  characterized  by 
high  sensitivity  to  initial  condition  errors.  For  a  small 
tolerance  level  (e  1),  the  model  prediction  skill  is 
represented  by  the  e-folding  time.  Both  the  mean  and 
variance  of  FPT  increase  with  increasing  intensity  of 
the  attractor  fluctuations.  The  mean  FPT  is  longer  in  the 
nonlinear  than  in  the  linear  regimes.  The  model  predic¬ 
tion  skill  should  weakly  depend  on  subgrid-scale  pa¬ 
rameterization  schemes;  however,  neglect  of  the  noise 
scale  (5  =  0)  leads  to  an  incorrect  asymptote  for  the 
mean  FPT  when  the  initial  error  is  infinitesimal.  For  a 
large  tolerance  level,  the  mean  FPT  is  no  longer  ex¬ 
ponentially  decaying  with  the  initial  error,  and  is  always 
greater  than  the  e-folding  time. 

The  high  predictability  (A  <  1)  is  characterized  by 
the  power-law  decay  of  the  mean  and  variance  of  FPT 
with  decay  exponents  (1  —  A)  and  2(1  —  A),  respec¬ 
tively,  if  z  »  1  and  A  ->  1 .  The  prediction  skill  strongly 
depends  on  the  noise  level  5  and  the  subgrid-scale  pa¬ 
rameterization  schemes.  Notice  that  both  the  mean  and 
variance  of  FPT  weakly  depend  on  the  values  of  infin¬ 
itesimal  initial  errors. 

A  nonlinear  perspective  of  forecast  error  is  important 
in  finding  an  exact  functional  correlation  between  the 
quantitative  model  skill  and  the  physical  factors.  In  gen¬ 
eral,  nonlinear  processes  limit  the  growth  of  forecast 
error  and  causes  a  wide  dispersion  of  individual  pre¬ 
dictions.  Averaging  over  the  ensemble  of  initial  uncer¬ 
tainties,  the  nonlinearity  can  either  increase  or  decrease 
the  mean  FPT. 

The  model  prediction  skill  strongly  depends  on  the 
position  of  the  predicted  attractor.  Low-  or  high-  (local) 
prediction  skill  should  be  locally  detected  through  PDP 
analysis  on  different  fragments  of  the  attractor  trajec¬ 
tory.  However,  the  real  model  skill  is  obtained  through 
averaging  over  the  local  prediction  skills  along  the  ref¬ 
erence  trajectory.  The  high  (low)  predictability  can  be 
interrupted  by  randomly  distributed  bursts  of  low  (high) 
predictability.  Such  a  temporal  intermittency  of  the  pre¬ 
diction  skill  should  be  investigated  in  more  complete 
model  than  the  Nicolis  model.  This  is  not  discussed  in 
this  paper. 
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APPENDIX  A 

Iterative  Process  for  Solving  (5.5)  and  (5.6) 


Taking  (5.6)  as  an  example,  the  iterative  process  starts 
from  the  decomposition  of  the  operator 


L- <*6. 

2  at;  o 


(A.l) 


into  two  parts 


L  =  L,  +  iiL2,  (A.2) 

where  /z  is  the  varying  parameter  ranged  from  0  to  1 . 
The  operator  L,  is  chosen  such  that  the  equation 


Lx  t=— 1  (A.3) 

A  1  r  '  - 

has  an  analytical  solution. 

Then,  we  introduce  an  additional  auxiliary  function 
W  =  dT/dfjt  satisfying  the  following  equation: 

LlW  +  /jlL2W  +  L2(t)  =  0.  (A.4) 

•  i ... 

The  following  iterative  process  is  used  to  solve  (A.4) 
(Ivanov  et  al.  1999): 

Lii r)l0i  ff  “I,  (A.5) 

<*>»  =:(T^'r>i  +  A/iWt'-H, 

/  =■  1,  ....  7  (A.6) 

+  L2(t)1'1  +  /*“iL2Wt'-»  =  0,  (A.7) 

=  Ajiu,  A  jj.  =  (A. 8) 


where  I  is  the  maximum  number  of  iterations.  The  it¬ 
eration  stops  when  /j3‘]  equals  1.  Here,  (t)[,]  is  the  ap¬ 
proximate  solution  of  (5.5)  with  accuracy  of 
The  variance  of  FPT  is  calculated  using  the  same  it¬ 
erative  procedure. 

Let  us  illustrate  the  iterative  procedure  (A.5)-(A.8) 
by  calculating  the  mean  FPT  for  A  =  1  and  1=3: 

t 
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(A.9) 

LjWm 

(A.  10) 

W[11 

=  -L2<t)>°L- 

(/x[1i  +  A  pl)L2W«», 

(A. 11) 

L,Wm 

=  -L2<T>“»  - 

ApL2Wm 

-  (yu,i2i  +  A/t^W™, 

(A.  12) 

L,W™ 

=  -^(r)™  - 

AptL2Wm  -  ApL2W™ 

-  (/xt3i  +  A pl)L2W™, 

(A.13) 

<T)foj  =-/?(-,  2  |, 

W[01  =  %3I, 


3  q 


W 


8  8 
e  & 


8  8 
e  £o 


W[11  =  3  +  —  B2A^R k,  T.  4  • 


9  q2 


e  £o 


S  5 


W[2)  =  e,  T.  3  +  7T-B2A ^R  4 


9q 2 


8  8 
e  £o 


+ 


45$2 


e  £o 


*3*!  g,  7-5  1. 


3  q2 


e  £0 


8  8 


W‘3>  =  — B.i?  3  +  —  B2A^R  4 


3?2 


5  8 


e  £0 


8  8 


36 

135<?2 


e  £0 


S’  S 


16 


«3ii;,!,5+-r:Mb7,6 


405g2 


S  £o 


S  8 


where 


R(.Zi,  z2,  l )  =  ln'^  -  ln'z2,  Bn  =  , 

n  =  1,  ....  4. 

APPENDIX  B 

Averaging  over  Initial  Uncertainty 

a.  Small  amplitude  error  and  small  tolerance  level 


Let 


tmin  —  t  t  tmax  =  t  _l_  t 

bO  bO  b>  bO  bO  '  b> 

0  =  lo  =  0.5(f5“  +  £r)- 

i  .  •-  -  %£i  ’  .U*£  5s  !>ii4 


is  given  by 


<T>  = 


1 


q2(y  -  1) 


In 


tma 

bO 


+  1 


+ 


tmin 

bO 


lnl 


where  Eqs.  (A.9)-(A.13)  use  the  reflecting  boundary  con¬ 
dition  at  £0  =8  and  the  absorbing  boundary  condition  at 
£0  =  e,  respectively.  The  solutions  are  written  by 

1  JE  & 


tmir 
bO 

t max  _  tmin  \  tma 

bo  bo  \bo 
Rewrite  (B.l)  and  (7.2)  into 


(B.l) 


«T(|0)»  = 


q\y  -  l) 


<T(lo)>  = 


q2(y  -  1) 


In 


In 


e 

£max 

SO 


£m£ 

SO 


+  1  + 


+  In 


1  -  0 


1  +  0 


ln(0) 

(B  2) 

.  (B.3) 


A  comparison  between  (B.3)  and  (B.4)  leads  to 
«T(|0)»  >  <r(|0)>  if  0  £  0.9. 

b.  Small  amplitude  error  and  large  tolerance  level 
For  this  case,  £0/s  <  l,a  constant  A, 

4 


A  = 


3<?2(A  -  1)|?2 


ge  gW  +  +  g4e4 

q 4  q 6  q% 


is  added  to  both  (B.2)  and  (B.3).  The  analysis  becomes 
the  same  as  the  previous  case  (i.e.,  the  small  amplitude 
error  and  small  tolerance  level). 


c.  Large  amplitude  error  and  large  tolerance  level 
Use  of  the  asymptotic 


S  G 

ln(  -f- 1  =  —  -  1  as  —  -4  1, 

(So/  So  So 

and  average  of  all  terms  of  (7.4)  over  the  ensemble  of 
initial  errors  leads  to 


<t(|0)>  = 


«t(|0)»  = 


1 


q\y  -  1) 


G  + 


4 g2e2A. 


:.j£uJ3uJ 


i 


q2(y  -  l) 


G  + 


4  g2e2A2 


If; 


2 

max  / 


1  -  7^(1  +  20  +  02) 

4  £2 


+ 


4g2e2A3 


3  q6 


1  -  -^(1  +  20  +  02) 

4e2 


3?4 


1  -  !%(1  +  0  +  02)  + 

3  e2  3 q6 


& 


3,1  ~  5(1  +  p2) 


+  . . . 


+  "t 

(B.4) 
,  (B-5) 


t 
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where 


Comparison  between  (B.4)  and  (B.5)  leads  to 
«T(|0)»  <  <T(|„)> 
for  any  value  of  f3. 
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